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I. SCOPE 

In a previous paper which we shall call Paper I, 
we showed that one can uniquely define the Lemaitre- 
Tolman cosmological model 2, 3] by specifying an initial 
density profile (i.e. the mass-density as a function of the 
radial coordinate at an initial instant t = t\) and a final 
density profile. The formulae defining the L-T functions 
E(M) and ts(M) (where M is the active gravitational 
mass, used here as a radial coordinate) are implicit but 
unique, and can be solved for E and ts numerically. (For 
definitions of E and t B see eqs. (12 . 1[) and (|2.2[) .) We 
also worked out a numerical example in which a galaxy- 
cluster-like final profile was created out of an initial pro- 
file whose density amplitude and linear size were small. 

In the present paper, we develop that study for new 
elements: we show that instead of a density distribution, 
one can specify a velocity distribution (strictly speak- 
ing, this is i? j t/M 1 / 3 — a measure of the velocity) at 
either the initial instant or the final instant or both. We 
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prove a theorem, analogous to the one proven in paper I: 
given the initial and the final profile, the L-T model that 
evolves between them is uniquely determined. We also 
show how to adapt the initial and final density profiles to 
the astrophysical data more precisely than it was done in 
paper I. We provide numerical examples of L-T evolution 
between an initial profile (of density or velocity) consis- 
tent with the implications of the CMB measurements, 
and a final profile that corresponds either to a galaxy 
cluster or to a void. 

The paper is arranged as follows. We recall the ba- 
sic properties of the L-T model in sec. 2. In sec. 3 we 
describe how the final profile can be adapted to obser- 
vational data more exactly than in Paper I. In sec. 4 
we find the implicit formulae to define the L-T functions 
E(M) and is(M) when the initial state is specified by 
a velocity distribution and the final state by a density 
distribution. In sec. 5 we do the same for both states 
being specified by velocity distributions. In sec. 6 we 
deduce the amplitude of the velocity perturbation at t\ 
allowed by observations of the CMB radiation. In sec. 
7 we specify our choices of density and velocity profiles 
for numerical investigations. Sees. 8 and 9 contain the 
presentation of numerical results in several figures. The 
results are summarized and interpreted in sec. 10, and 
brief conclusions follow in sec. 11. 

The approach presented in this and in paper I is more 
suited to astrophysical practice than the traditional ap- 
proach, in which we first try to deduce the initial state 
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of the Universe from various kinds of data, and then pro- 
ceed to calculate the evolution of the cosmological model, 
in order to compare its final state with the observations 
of the current state of the actual Universe. Our approach 
allows one to make simultaneous use of the data on the 
initial and on the final state of the Universe — the real 
astronomical data is indeed such a mixture. 

Naturally the L-T model can describe the actual as- 
trophysical process of structure formation only approx- 
imately. The obvious limitation is spherical symmetry, 
in consequence of which we cannot take into account the 
rotation of the objects formed. Thus, no matter how well 
we reproduce the profile of a galaxy or cluster, the L-T 
model will continue to evolve, and the profile will look 
quite different after 10 6 years or less, whereas real galax- 
ies and clusters are fairly stable over 10 9 years. However, 
we hope that the method presented here will be a starting 
point for generalisations that will be done once (and if) 
more general exact cosmological models are discovered. 



II. BASIC PROPERTIES OF THE 
LEMAITRE-TOLMAN MODEL. 

The Lemaitre-Tolman (L-T) model @, S] is a spheri- 
cally symmetric nonstatic solution of the Einstein equa- 
tions with a dust source. Its metric is: 



Parabolic case, E = 0: 



ds 2 = dt 2 - 



R, r ' 



1 + 2E(r 



■dr 2 - R 2 (t,r){d$ 2 + sin 2 i?d<p 2 ), 

(2.1) 

where E(r) is an arbitrary function (arising as an integra- 
tion constant from the Einstein equations), i?, r = dR/dr, 
and R obeys 



R/ = 2E + 2M/R + lAR 



(2.2) 



where A is the cosmological constant. Eq. (|2.2p is a first 
integral of the Einstein equations, and M(r) is another 
arbitrary function that arises as an integration constant. 
The mass-density is: 



up 



2M, r 
R 2 R, 



where k 



IttG 



(2.3) 



See [J for an extensive list of properties and other work 
on this model. In the following, we will assume A = 0. 
Then eq. (|2.2|) can be solved explicitly. The solutions are 
the following. 



Elliptic case[2J, E < 



R(t,r) 



M 
2E 



(1 — COST?), 



{-2Ef/ 2 
r\ - sin?? = — (t - t B (r)). 



M 



(2.4a) 



(2.4b) 



where r\ is a parameter; 



R{t,r) 
Hyperbolic case, E > 



-M(t-t B (r)) 2 



1/3 



M 



(2£;)( 3 / 2 ) 
sinhr? - r\ = — (t - t B (r)), 



(2.5) 



(2.6a) 



(2.6b) 



where t B (r) is one more arbitrary function (the bang 
time) . Note that all the formulae given so far are covari- 
ant under arbitrary coordinate transformations r = g(r'), 
and so r can be chosen at will. This means one of the 
three functions E(r), M(r) and t B (r) can be fixed at our 
convenience by the appropriate choice of g. 

The Friedmann models are contained in the Lemaitre- 
Tolman class as the limit: 



t B = const, 



E\ 3/2 /M = const, 



(2.7) 



and one of the standard radial coordinates for the Fried- 
mann model results if, in addition, the coordinates in 
- ([23]) are chosen so that: 



M = M r 6 , 



(2.8) 



where Mq is an arbitrary constant. This implies E/r 2 — 
const := — k/2, k being the Robertson- Walker curvature 
index. 

It will be convenient in most of what follows to use 
M(r) as the radial coordinate (i.e. r' = M(r)) because, 
in the structure formation context, one does not expect 
any "necks" or "bellies" where M, r = 0, and so M(r) 
should be a strictly increasing function in the whole re- 
gion under consideration. Then: 



Kp = 2/(R 2 R, 



M 



6/(i? 3 ). 



M ■ 



from which we find 



R 3 (M) - Rl 



3 

47T 



M 



Mo 



du 

p(u)' 



(2.9) 



(2.10) 



where Rq and Mq will commonly be zero. 

In the present paper we will apply the L-T model to 
problems of a similar kind to that considered in Paper 
I: Connecting an initial state of the Universe (defined 
by either a mass-density distribution or a velocity dis- 
tribution) to a final state (also defined by one of these 
distributions) by an L-T evolution, and in particular to 
the formation of galaxy-cluster-like and void-like objects 
out of initial perturbations of density or velocity that 
are small in amplitude and in some cases small in mass 
compared to the final object. 
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III. MODELLING THE FINAL DENSITY 
PROFILE. 

We will now incorporate in our models the observa- 
tional data on mass distribution in galaxy clusters in a 
more detailed way than in Paper I. 

The quantities of interest in the profile are the follow- 
ing: the maximal density (with the shapes we assume 
below and later, this will be at the center of the object), 
the radius of the object (assumed spherically symmetric), 
the mass of the object, the average density of the cosmic 
background, and the compensation radius (defined be- 
low). 

We define the following parameters. 

• M m — the mass of a galaxy cluster, out to some 
radius; to be taken from astronomical tables. 

• R m — the radius within which the mass M m is con- 
tained; also to be taken from astronomical tables. 
In fact there may be two or more M m -R m pairs 
available for some clusters (see section [VII C ip . 

• p(M m ) — pbd — this is the geometrical definition of 
the mass M m (M — M m at that value of R = R m 
at which p equals a certain specified multiple of the 
background density). [25[ We do not assume any 
value of d at this point yet, except that 1 < d < 

fmax IPb- 

• The compensation radius R = R c , at which the 
total mass within R = R c is the same as the back- 
ground mass would be if no inhomogeneity were 
created. This is needed to let us know roughly 
where our inhomogeneity is matched to the Fried- 
mann background. See subsection IIIIAI 



A. Compensation Radius &: Mass 

The matching of contained mass at R c is not suffi- 
cient for a Swiss-cheese type matching, as we have not 
required the time since the bang to match up at this loca- 
tion, and therefore is technically the wrong definition [2t|. 
However, it has the advantage that it can be calculated 
knowing only the density profile at a given time, and 
under the circumstances used here is a fairly good esti- 
mate. All our models are then calculated out to M c , i.e 
< M < M c , 

We can put upper and lower limits on the compensa- 
tion radius and mass. For a condensation of measured 
mass M m , the M c value obviously cannot be less than 
this. In fact the region around the visible condensation, 
though of low density is of large volume, and will add no- 
ticeably to the mass to be 'compensated for'. Therefore, 



M r 



> M„ 



Rc 



> 



( 3M„ 



\47rp, 



6, to day 



1/3 



(3.1) 



where eq. (|7.2[) gives p^today = 4.075 x 10 30 g/cc. in 
the chosen parabolic background, or more conveniently 



Rc. 



> 1.58 x 10- 4 {M rn /M Q ) 1/3 in Mpc. (3.2) 



On the other hand, the observed average separation 
of condensations puts an upper limit on R c . Since the 
contents of the universe are a mixed bag of galaxy clus- 
ters, rich clusters, superclusters, field galaxies, voids, 
walls, etc, the average separation of rich clusters say is 
not meaningful for this purpose. So we instead argue 
that there are around 3 x 10 6 large galaxies & 3 x 10 7 
dwarf galaxies within 10 9 light years. At say 10 11 M Q 
and 10 10 Mq, respectively, this gives a mean density of 
around 4.97 x 10 9 M /Mpc 3 = 3.36 x 10" 28 kg/m 3 = 
0.0826pb ! today Therefore an object of mass M m in M 
should on average occupy a volume of radius about 



i? c , max = 3.64 x 10" 



'(MjMg) 1 ' 3 in Mpc. (3.3) 



For a galaxy of 10 11 Mq, these two limits are 0.735 < 
R c < 1.69 Mpc, whereas for an Abell Cluster of 10 15 M©, 
these two limits are 15.8 < R c < 36.4 Mpc. 

For a void, the interior density will not be zero, but is 
not well known, and the radius is of the order of 60 Mpc. 
Only by including some of the galaxies in the surrounding 
walls can one bring the average density up to the back- 
ground value. So R c > 60 Mpc, which at background 
density gives M c > 5.5 x 10 16 M . 



B. An example of fitting a profile. 

Since little is known about mass-distribution within 
galaxy clusters, we cannot attempt to model any actual 
cluster to any significant accuracy. However, we wish to 
show that such modelling is in principle possible. There- 
fore, for the beginning, we will use the profile from Pa- 
per I and will show below that its free parameters can be 
made equal to some observed/observable quantities. 

Let pb be the background density. We choose the pro- 
file at t 2 to be 



p(M) = B 2Pb e-^^ 2 



(3.4) 



where B 2 and p 2 are free parameters to be adapted to 
the observational constraints. 

Given two pairs of (M m ,R m ) data, (M a ,R a ) & 
(Mb, Rb), we have 



R 3 



3^2 



-erfi(M//z 2 ) 



(3.5) 



8^ B 2 pb 

for each of them, and each can be solved for B 2 , so that: 

3/i 2 erfi(M a //i 2 ) _ _ 3p 2 crft(M b / p 2 ) 

8V^ RlPb 2 ~ 8V¥ Rlpb ' [ ' ' 

This is solved numerically for p: 2 , and B 2 follows. 
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The compensation radius R c is determined by the con- 
dition that the mean density out to M equals the back- 
ground density: 



Pb = 



2B 2 p b M c 



4ttR(M c ) 3 /3 ^ p 2 crti{M c /p 2 ) ' 



(3.7) 



This is solved for M c numerically and R c follows from 



C. Profiles inspired by astronomy. 

We looked for density profiles that are considered real- 
istic by astronomers. As it turned out, there is no general 
agreement as to which profile best describes observations, 
and no generally accepted definition of the radius of a 
galaxy cluster exists. However (see the Acknowledge- 
ments) the following 'universal profile' is one of the more 
commonly used formulas for density vs. distance profiles 



p(0 = PbS 



1 



(C/C s )(i + C/C) 5 



(3.8) 



where pb is the average density in the Universe, S is a 
dimensionless factor and £ s is a scale distance. (This is a 
Newtonian formula, in which £ is the Euclidean distance.) 
According to the authors of Refs. d, [UGH, this profile 
applies for £ changing by two orders of magnitude. 

For our procedure, we need the density given as a func- 
tion of mass, which cannot be found analytically for the 
above profile. [27] The calculation p(Q -> M(() -> 
C(M) — > p(M) can always be done numerically, but it 
is more instructive to have exact explicit formulae. 

We therefore approximate M(() separately in the 
ranges ( <C ( s and C > Cs • 

For £ <C ( s a Taylor series around £ = up to £ 2 gives 



M(C) « 2irp b 5C s ( 2 := /3 2 ( 2 -> C « VM/0. (3.9) 
which converts the profile f|3.8[) to: 

1 



Pl {M) = Pb j35( s 



(3.10) 



Like the original profile, this one has the unpleasant prop- 
erty that the density becomes infinite at M — > 0, so we 
modify it to 



p 2 (M) 



Pbf38( s 



(e + VM) (l + ^g)' 



(3.11) 



where e is small compared to /3<5£ s , and this can be inte- 
grated to give R(M), 

When ( ^S> ( s , the approximation 



M(C) » lirpbSCs In (1 + C/C-) ~ 7 In (1 + C/G) , 



C = C s (c Mh i) , 



(3.12) 



substituted into (13 . 8[) gives 



p4,(M) = p b - 



Mh - 1 + , 



,-2M/ 7 



(3.13) 



where v has again been added into the denominator to 
permit a non-divergent central density. The correspond- 
ing R 3 (M) is still elementary. [H 

The profiles (|3.1ip and (|3 . 1 3|1 were the starting point 
for our considerations. The profiles actually used in nu- 
merical examples were various modifications of these (see 
Sec. IVIip . done in order to better fit our profiles to ob- 
servational data. 



D. The initial density profile. 

We assume that the condensed object (model of a 
galaxy cluster) was created out of a small localized ini- 
tial density or velocity perturbation, superimposed on 
the homogeneous spatially flat Friedmann background. 
During the evolution, the perturbation was increasing in 
density amplitude and in mass, and thus was effectively 
drawing more mass into the condensation region out of 
the surrounding homogeneous region[29l|[30l|. 

At t = ti, the profile need not be compensated. For 
an uncompensated profile, it is good if it is localized (i.e. 
the perturbation is zero for M > Mi, where Mi is the 
assumed mass of the initial perturbation). Then the def- 
initions of the radius and mass of the perturbation are 
straightforward. 

We choose for this profile the cosine shape 



I M\ - / + cos(ttM/Mi)) + p b for M e [0, Mi]; 
P[ >-\pb = const for M > Mi 

(3.14) 

with Bi > 0. Then p(Q) = pb(l + A) is the maximal value 
of density, with Ap/pb = A = 2Bi/pb being the density 
amplitude from astronomical data, while Mi is the total 
mass of the initial perturbation. 

The radius R(M) is given from eq. (|2.10p by 



i? 3 (M) 



x arctan 



3Mi 



2ttVp6 (Pb + 2B 1 ) 

( 7rM 



Pb 



tan 



p b + 2B 1 \2Mi 



for M < Mi, 

(3.15) 



R^M) = 3Ml + (M - Mi) 

Any/pt (p b + 2Bi) 4ttP6 (3.16) 

for M > Mi, 



where the first term on the right in (|3.16[) is the common 
value of R 3 at M = Ml. 

We did one run with this kind of initial profile (pi4 
in Table IVip in order to demonstrate more clearly the 
point we made in Paper I: that the mass of the "seed" 
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structure existing at time t\ can be much smaller than 
the mass of the galaxy cluster into which it will evolve. 
The other initial profiles are either flat (see Table IVIII 
and Sec. [X]fbr reasons why) or go into the background 
only asymptotically. 



IV. THE EVOLUTION FROM A GIVEN 
VELOCITY DISTRIBUTION TO A GIVEN 
DENSITY DISTRIBUTION. 

The quantity that is a measure of the velocity distri- 
bution of the dust in an L-T model is 



Rt/M 



1/3 



(4.1) 



This is constant in any Robertson- Walker model, so its 
nonconstancy is a measure of the velocity inhomogeneity. 

Suppose we wish to adapt an L-T model to a given 
initial velocity distribution b = b\{r) at t = t\, and to 
a given density distribution p = P2(f) at t = fa- This 
is a different set of data from the one considered in pa- 
per I, and so the existence of such an evolution has to 
be proven. The functions appearing along the way are 
different, but the overall mathematical scheme is essen- 
tially the same. As before, we will mostly be using the 
mass M as the radial coordinate, and in each case we will 
assume that at the initial instant t\ the configuration is 
expanding, so 



R t {t u M) > b 1 (M) > 0. 



(4.2) 



Analogous reasoning can be carried out for matter that 
is initially collapsing, but this situation is just covered 
by the time reverse of the method given here and is not 
relevant for the problem of structure formation. 



A. Hyperbolic evolution. 



We have, for E > 0: 
M 



R,t = 



'2ER 



V / (1 + 2ER/M) 2 - 1. 



If 



R,t/M 



1/3 



(4.3) 



(4.4) 



is the data, then denoting R\t=t t = Ri and solving the 
above for fii we find 



Ri = 



2M 1 / 3 



h 2 -2E/M 2 / 3 



(4.5) 



and so from the evolution equations 

rji = arcosh(l + 2ER X /M) 

'b 1 2 + 2E/M 2 / 3 



= arcosh 



6i 2 - 2E/M 2 / 3 



(4.6) 



\ 



b x 2 + 2E/M 2 / 3 
b x 2 - 2E/M 2 / 3 



1 — arcosh 



[2E) 3 ' 2 
M 



(h - t B ) 



h 2 + 2E/M 2 / 3 \ 
b x 2 -2E/M 2 / 3 ) 



(4.7) 



This is the equation that will be used as initial data at 
t = t t . 

The equation at t = t 2 will be the same as eq. (3.3) in 
paper I for i = 2 viz: 



(1 + 2ER 2 /M) - 1 - arcosh(l + 2ER 2 /M) 

3 / 2 (4-8) 



(2EY 
~M 



(*a - t B ) 



Here, the quantity R 2 = R(t 2 ,r) is calculated from the 
given p 2 (r) — p(t 2 , r) using eq. (|2. 10|) . 

Just as in paper I, we introduce the symbols 



x = 2E/M 2/3 , a 2 = R{t 2l M)/M 1/3 . 



(4.9) 



The set of equations to determine E(M) and ts(M) is 
then (from $t7$ and iji^l): 




— 1 — arcosh 



= x 3/2 (h - t B ) , 



(4.10) 



Y(l + a 2 x) 2 — 1 — arcosh (1 + a 2 x) — x 3 ^ 2 (t 2 — ts) , 

(4.11) 

From (14. lip we now find tB'- 



t-B — fa 



r 3/2 



1 + a 2 x) — 1 — arcosh(l + 022;) 



(4.12) 

and substituting this in (|4.10l) we obtain the following 
equation to determine x: 



$h{x) = 0, 



(4.13) 



where 



$jr(s) :=^{l + a 2 x) -1- 




, , , ,'bi +x 

arcosh(l + a 2 x) + arcosh — = 

M -x 



r 3/2 



[fa - h) • 



(4.14) 



We will use the functions b% (M) and a 2 (M) implied 
by the assumed R. t i(M) and p 2 (M) to find E(M) and 
ts(M), and then to find p(ti, M). This will tell us about 
the relative importance of the velocity and density dis- 
tributions for structure formation. In particular, with b\ 
= const, the p(t\,M) will tell how big the initial density 
inhomogeneity has to be when the initial velocity distri- 
bution is exactly homogeneous, while the final structure 
is given. 
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1. Conditions for solutions to exist 



Eq. (|4.19p may be rewritten as 



Now we have to verify whether the equation (x) = 
determines a value of x. We see that &h(x) is defined in 
the range 



< x < bi 



(4.15) 



It is easily seen that (0) = 0. For x — > b x 2 , the second 
and fourth term in (|4. 14[) tend to (— oo) and to (+oo), 
respectively, while all the other terms are finite. Since 



lim 

c— >b\- 



arcosh [{b x 2 +x) / {b x 2 - x)] 
^/[(b 1 2 + x)/(b 1 2 -x)] 2 -l 

the second term is dominant. Hence: 



lim $jj(x) = — oo. 



= 0, (4.16) 



(4.17) 



We also have 



a 2 



.3/2 



2h 



\/2 + a 2 x (b x 2 -x)' 
\fx\ H {x). 



2^ 



(4.18) 



Since X H , X = -±a 2 5 / 2 /(2 + a 2 xf' 2 - ^/{b^ - x) 3 < 0, 
the function Xh(x) is strictly decreasing and may have 
at most one zero. Hence, $>h,x can be positive anywhere 
only if Aff(0) > 0, i.e. if 



*2 



tl < 



V2 



-a 2 



3/2 



3V 



(4.19) 



Since we consider the evolution from t\ to t 2 , it is natural 
to assume t 2 > t\. Then, (|4.19[) can be fulfilled only if 



2/a 2 < bi 



(4.20) 



This inequality is easy to understand if we use (|4.5[) for 
b\ 2 (with R = i?i), and recall the definition of a 2l eq. 
(TO]) . Then P~2T)]) is equivalent to 



2 2 
R 2 Ri 



2E 



(4.21) 



This in turn implies that R 2 > R\\ otherwise there will 
be no E > evolution from R\ to R 2 . Thus, in passing, 
we have proven the intuitively obvious: 

Conclusion: If matter is expanding at t = t\, then an 
E > L-T evolution from R{t u M) to R(t 2 ,M) will exist 
only if R(t 2 ,M) > R(t 1} M)n 

Note that R(t 2 ,M) > R(t x ,M) does not imply 
p(t 2 ,M) < p(ti,M), contrary to the intuitive expecta- 
tion, and in contrast to the Friedmann models. This 
point is discussed in Appendix [Bl 



a 2 > 



9 / 4 



1/3 



which means that between t\ and t 2 the model must have 
expanded by more than the E = model would have 
done. (Applying the definitions P~4"]) and ([H?]) to (1231) 



and (|4~5| for the E = case shows ti - 4/(3&i ) = t B 
and a 2 = [(9/2(t 2 -ts) 2 ] 1 / 3 .) 

If (|4.19|) is fulfilled, then &h(x) starts as an increasing 
function in some neighbourhood of x = 0, beyond which 
there is exactly one maximum, and finally it decreases to 
— oo at x — b\ 2 , thus passing (only once) through a zero 
value somewhere between the maximum and x = b\ 2 . 
Hence, eqs. (I4.19|) and (|4.20p are a necessary and suffi- 
cient condition for the existence of an E > L-T evolu- 
tion connecting the initial state (defined by the velocity 
distribution (|4.4p ) to the final state defined by the den- 
sity distribution p 2 (M). As before, the resulting model 
will have to be checked for the possible existence of shell- 
crossings and regular maxima or minima. 



B. Elliptic evolution; the final state in the 
expansion phase. 

Assuming that the initial state at t = t± is in the ex- 
pansion phase of the Universe, we find, similarly as in 
(|43]) . that with E < 



R.t - 



M 



^2ER 



y/1 - (1 + 2ER/M) 2 , 



(4.22) 



and then, using (|4.4[) we obtain (|4.5p once again. Since 
this time E < 0, the inequality 2M^ 3 /R + 2E/M 2 / 3 > 
has to be fulfilled, but this is identical to 2M/R + 2E > 
that has to hold for an E < evolution anyway. In place 
of (|4"l)j) - (|4~7)) we obtain at t = t x 



r\ = arccos 



2/3 



b 2 + 2E/M 
bi 2 - 2E/M 2 / 3 



(4.23) 



arccos 



h 2 + 2E/M 



f2/3 



h 2 - 2E/M 2 / 3 

(-m 3/2 ,. , , 
(h-t B ). 



\ 



i - 



b x 2 + 2E/M 2 / 3 
b Y 2 - 2E/M 2 / 3 



(4.24) 



Since we assumed the final state at t = t 2 to be in the 
expansion phase, too, the equation at t = t% is the same 
as eq. (3.15) in paper I for i = 2: 



arccos (1 + 2ER 2 /M) - yfl - (1 + 2ER 2 /M) 2 



{-2Ef/ 2 



M 



(*2 - t B ) 



(4.25) 
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We introduce a 2 by (|4.9[) and 

x = -2E/M 2/3 , 
then (1424)) and (|4~25"i) become 



h 2 -x 
b 2 + x 



7 2 x 2 

61 — a; 



6i 2 + .T 



(4.26) 



z 3/2 (ti - t B ) , 



(4.27) 



arccos (1 — 022;) — yl — (1 — a 2 a;) 2 = x 3 ^ 2 (£2 — is) ■ 

(4.28) 

From (|4.28l) we have 



tB — t 2 



c 3/2 



arccos (1 — a2x) — y 1 — (1 — g^x 



(4.29) 



and substituting this in (|4.27[) we obtain 

- 0, (4.30) 

where 



arccos (1 — 02a;) — arccos ■ „ 



b x 2 - x 



bi +x 



x 3 / 2 (h-h). 



(4.31) 



The terms containing b\ do not impose any extra restric- 
tion on x apart from x > that follows from the defini- 
tion (I4.26p . Since 02 > and x > 0, the restriction on x 
imposed by the other terms is the same as in paper I: 



x < 2/a 2 , 



(4.32) 



which follows from 2M/R + 2E > (x = 2/a 2 corre- 
sponds to the given M-shell being exactly at the maxi- 
mum of expansion at t = ti). 



1. Conditions for solutions to exist 

The calculations are similar to those in the previous 
subsection, they are described in Appendix IC II Eq. 
(I4.30P has a solution if and only if the following two in- 
equalities are obeyed: 



V2 o/n 4 

2_ 1>_ T a2 36? : 



(4.33) 



which is the opposite of (|4.19[) . and 
ti - *i < 



(oa/2) 



3/2 



6iV2a 2 
a 2 ^i 2 /2 + l 



arccos 



/ q 2 fei 2 /2-l 
Va2&i 2 /2 + l 



(4.34) 



The solution is then unique. The second inequality is 
equivalent to eq. (3.22) in paper I, as can be verified by 
writing (|4.5[) as b\ = \j2ja\ — x and using x = 2ja%. 

It follows that the right-hand side of (|4.34|) is always 
greater than the right-hand side of (|4.33|) . as proven in 
Appendix A to paper I, so the two inequalities (|4.34|) and 
(I4.33[) are consistent. It follows further, as shown in pa- 
per I, that the equality in (|4. 34[) means that the final state 
at t-z is exactly at the maximum expansion. Note that, in 
consequence of b\ > 0, ai > and arccos (anything) < n, 
the right-hand side of (|4.34[) is always positive, and so, 
given ti, values of ti obeying (|4.34|) always exist. 

The set of inequalities (|4.33[) - J4.34[) constitutes the 
necessary and sufficient condition for the existence of 
an E < L-T evolution between the initial state at 
t\ specified by the velocity distribution (|4.4[) and the 
final state at ti specified by the density distribution 
p(ti,M) = pi(M), such that the final state is still in 
the expansion phase. 



C. Elliptic evolution; the final state in the 
recollapse phase. 

We assume now that the initial state at t = t\ is again 
in the expansion phase, so eqs. (|4.22|) - (|4.24|) still apply. 
However, the final state at t = ti is now assumed to be 
in the recollapse phase. Instead of eq. (|4.25[) we then get 



7T + arccos (-1 - 2ER 2 /M) + \Jl-{l + 2ER 2 /MY 
_ {-2Ef 2 



M 



{ti - t B ) ■ 



(4.35) 



Introducing again a 2 by (I4.9|) and x by (|4.26[) . we ob- 
tain the set of equations consisting of (|4.27p and 



7r + arccos (—1 — a 2 x) + yl — (1 — aixY 
= x 3/2 (ti-t B ). 
The solution for is is now 
tB = ti 
1 

~ Z3J1 



(4.36) 



7r + arccos (— 1 + a 2 x) + y 1 — (1 — a 2 x 



and substituting this in (|4.27[) we get 

&c(x) = 0, 

where 



(4.37) 
(4.38) 



bi- 



bi 



+ (1 - a 2 x) 



(hi -x\ 

tt — arccos — 5 + arccos ( — 1 + a 2 x) 

\bi+xj 

,3/2 



(4.39) 
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As before, we have x > by definition and x < 2/a 2 . 



of (|4.11|) we now have 



1. Condition for solutions to exist 

The details of the calculation are shown in Appendix 
C2l Eq. (|4.38[) has a solution if and only if 



h-U> (a 2 /2) 3/2 



0261/2 + I 



/ a 20l 2 /2-l 
Va 2 6i 2 /2 + l 



(4.40) 



and then the solution is unique. This is the necessary 
and sufficient condition for the existence of an E < 
L-T evolution between the initial state at t\ specified by 
the velocity distribution (|4.4[) and the final state at ti 
specified by the density distribution p(t 2 ,M) = p 2 (M), 
such that the final state is in the recollapse phase. 



D. A note on the meaning of parameters. 

It should be noted that the values of the time coordi- 
nate ti and t 2 , at which the initial and final states are 
specified, do not play individual roles in the calculations. 
The meaningful parameter is fa — ii); it is this value 
that determines E(M), and then the corresponding ts is 
calculated. The "age of the Universe" at the two instants 
then follows, being fa — ts) and fa — £b)j respectively. 
This conclusion applies as well to paper I; and the same 
is still true in the Friedmann limit, where x, a\, a 2 and 61 
are no longer functions of M , but just constants. Hence, 
the physical input data for the procedure, both here and 
in paper I, are the initial and final distributions of density 
or velocity and the time- interval between them, fa — ii)- 
The individual values of t\ and t 2 do not have a physical 
meaning. 




— 1 — arcosh 



b 2 2 



= x 3/2 fa - t B ) , 



b 2 = R t /M 



2/3 



(5.1) 



(5.2) 



Since the hyperbolic expansion cannot be reversed to be- 
come collapse, and since we assumed bi > (expansion 
rather than collapse at ti), it follows that 6 2 > is a 
necessary condition here. 

From (|4.10p and (|5.1| we obtain two equations 




1 — arcosh 



so x is the solution of 



Xh = 0, 



(5.3) 
(5.4) 





(b 2 2 + x\ 








W-x) 





r 3 / 2 (U - 



(* 2 - h) . 



(5.5) 



Note that (|4.5[) and its analogue at t 2 imply bi — x = 
2/di > 0, i = 1,2, i.e. x < bi 2 at both i = 1 and 
i = 2. This, together with x > 0, guarantees that the 
argument of arcosh will be always > 1, as it should be. 
Simultaneously, it guarantees that both square roots will 
be real. 



V. THE EVOLUTION FROM A GIVEN 
VELOCITY DISTRIBUTION TO ANOTHER 
VELOCITY DISTRIBUTION. 

For completeness, we shall also consider the case when 
both the initial and the final state are defined by veloc- 
ity distributions, even though this does not seem to be 
directly useful for astrophysics. 



A. Hyperbolic evolution. 

Eq. (TjnTD stiU applies. Defining x = 2E/M 2 / 3 as in the 
previous cases, we again obtain eq. (|4.10p . while instead 



1. Conditions for solutions to exist 

Details of the calculations are shown in Appendix IC 21 
For further analysis, we will need to investigate the 
derivative of Xh{x), it is 



Ax 



Xh{x) 



2b 2 



2bi 



(b 2 2 -x) Z (b, 2 - 



\/x^ H (3 



(t 2 - h) 



(5.6) 



We also have to know which of b\ , b 2 is greater. As 
expected, it turns out that in consequence of the assump- 



tions made (t 2 > t\ and i?., 



it=ti 



> 0, i.e. h > 0), bi 
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must be greater, or else the E > evolution between the 
two states will not exist — see Appendix [D] 
Therefore we take it for granted that 



h > b 2 , 



(5.7) 



i.e. that at every point the velocity of expansion at t 2 
must be smaller than at t\ at the same comoving coordi- 
nate M. With b\ > b 2 , we have 

x < b 2 2 , lim ^h{x) = +oo, 



x— »6 2 



(5.8) 



and, in a similar way as in Appendix [Dl it can be proven 
that now n' H (x) > for all x > 0. Consequently, (1h(x) 
may have at most one zero, and it will have one only if 
/if/(0) < 0, i.e. if 



1 



1 



t2 - h> -A^~^)- (5 - 9) 

Unlike in those cases where the profile was specified by 
density, the timc-diffcrcncc between the initial and the 
final state must be sufficiently large. The inequality be- 
comes more intelligible when it is rewritten as 

which means that an E > evolution between the two 
states will exist provided the velocity of expansion at t 2 
is greater than the velocity of expansion of the E = 
model, for which b 3 = 4/[3(t-t B )], giving an equality in 
(f5T0l) . 

Since the range of x is limited to (0, b 2 2 ), finding x 
numerically will pose no problems. 



B. Elliptic evolution, the final state in the 
expansion phase. 

Expansion at t 2 is equivalent to b 2 > 0. The set 
of equations to define x = —2E/M 2 / 3 now consists of 
(I4.27[) and its analogue at t 2 , obtained by the substitu- 
tion (&i,<i) — > (b 2 ,t 2 ). Hence: 



ts =t 



1 



3-3/2 
1=1,2, 

and x satisfies 

where 

Xx (x) ■= - 



arccos 



1- 



Xx(x) = 0, 



(5.11) 
(5.12) 




This time there is no additional limitation imposed on 
x by eqs. (|5.11[) - (|5 . 13[) — the square roots and the 
arccos will be well-defined for every value of x > 0. The 
relation between x and 6, is x = 2M 1 / 3 / (b 2 + l)] , 
and so indeed any value of x > is allowed, since 6; can 
have values between at maximum expansion and oo at 
the Big Bang, while Ri can be arbitrarily large or small, 
depending on the initial data. However, as always in the 
elliptic case, E > —1/2 has to hold, E — —1/2 indicating 
a maximum in the spatial section, and this will have to 
be checked while solving (|5.13[) for x. 

Finding x numerically will be no problem also here, 
since the natural variables to work with are j/j = 
(pi 2 — x) I {pi 2 + x) , both of which have the limited 
range (—1,1) for x > 0. 



1. Conditions for solutions to exist 

We have: 

Xx(0) = 0, lim xx (x) = -oo, (5.14) 



X'x(x) = \/x 



2b 2 



2bi 



Y(b2 2 - 

y/xfj,x(x). 



2 (*2 " h) 



(5.15) 



Xx(0) = 0, lim XxW = -oo. (5.16) 



It is known from the evolution equations (|2.4[) that if 
b\ > and t 2 > t±, then the velocity of expansion at t 2 
must be smaller than b\. However, just as in the previous 
cases, this can be proven also from the properties of the 
function xx (x) , and the proof is given in Appendix [E] 
Therefore we will take it for granted that 



b 2 < h 
in the following. We find 

Ab 2 4&i 



Hx(x) 



(b 2 2 +x) 3 (b^+x) 3 ' 



(5.17) 



(5.18) 



and so n' x (0) < 0, i.e. ^x(x) is decreasing in a neigh- 
bourhood of x = 0. It keeps decreasing as long as p! x is 
negative, i.e. as long as x < x m , where x m is the (unique) 
solution of n'x{x) = 0: 

. 5/3 . 5/3 
/, , n1/3 01 —02 , a/3 

x m = (bib 2 ) — -- — -— r = (hb 2 ) 1 x 



(5.13) 



(bi 4/3 + 6i6j 



6 1 1 / 3 -6 2 1 / 3 



1 / 3 + 6 1 2 / 3 6 2 2/3 + i 1/3 a 2 + 6 2 4/3 



(5.19) 
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From here we find 



.1/3. 



t 2 i 2 
01 - 02 



(5.20) 



and so the value of (ixfam) is 



Mx(xm) =2 



/ 3 - 6 2 1/3 



(oi 2 - 6 2 2 ) 2 



6 2 4/3 -6i 4/3 



(5.21) 



| (ta - *i) < -| (*a - *i) < 



(in consequence of 62 < 61). For x > x m , fJ.' x (x) is pos- 
itive, and so Hx{x) is increasing, but lim^^oo Hx{x) = 
— 1(^2 — ti) is still negative. Hence, in summary, nx(x) 
is decreasing from the value 2/6 2 3 — 2/oi 3 — | (i 2 — ti) at 
x = (which may be positive or negative) to a minimum 
at x = x m which is negative, and then keeps increasing, 
but remains negative for all x > x m . Thus, from (|5.15[) 
it follows that the function \x is necessarily decreasing 
for x > x p > 0, where x p is the zero of x'x( x )- Con- 
sequently, the equation Xx{%) = will have a positive 
solution only if x p > 0, i.e. if there is a neighbourhood of 
x = in which xx(x) is increasing (recall: xx(0) = 0). 
Thus, the necessary and sufficient condition for the exis- 
tence of a positive solution of xx(x) = is ^x(0) > 0, 
i.e. the opposite of (|5.9p : 



(5.22) 



4/1 1 

which now implies that the expansion between t\ and £2 
must have been slower than it would be in an E = 
model. 



C. Elliptic evolution, the final state in the collapse 
phase. 

Finally, we shall consider the E < evolution between 
two velocity profiles for the case when the final state is 
in the collapse phase, b 2 < 0. 

For the initial state, we can still use eqs. (|4.3p - (|4.7|) . 
The only difference is that this time b 2 comes out nega- 
tive: 



02 = -1 



'2MV3 2E 



(5.23) 



R 2 M 2 / 3 ' 
but the analogue of (|4.5p has the same form as before: 

2M 1 / 3 



Ro = 



2E/M 2 / 3 



(5.24) 



Hence, 
ts =fa 



,3/2 



7r -(- arccos 



x + 6 2 



1 - 



(5.25) 



with (|5.11[) still applicable at ti, and from (|5.25[) . the 
equation to determine x = —2E/M 2 / 3 is: 



Xc(x) = 0, 



(5.26) 



where 



Xc{x) :=\ 1 




h x l - x 



hi 2 +x 



x 3 / 2 (h-h) 



Conditions for solutions to exist 



(5.27) 



We have 



Xc(0) = 2tt, 



lim xc(x) = -00, 



(5.28) 



and so xc( x ) = is guaranteed to have a solution for 
x > 0. The surprising result is thus: 

// the final state is in the collapse phase, so that 6 2 < 
then an E < evolution exists between any pair of states 
(tl,&l) and (t 2 ,b 2 ). 

Consistently with (|5.28|) . we obtain: 



Xc( x ) = 



2b ly /x 



2o 2 3 



-Vx(t 2 -<i), 



{h 2 + xY (b 2 2 +x)-^ 

(5.29) 

which, in consequence of b 2 < 0, is negative for all x > 0, 
and so the solution of xc(x) — is unique. 

The initial point for solving xc{x) = numerically 
can be conveniently found as follows: each square root 
and each arccos-term in (|5.27[) is not larger than 1 at all 
values of x > 0, and where two of them equal 1, the two 
others are strictly smaller than 1. Hence, their sum is 
strictly smaller than 4, and 



X c(x) <4 + ir-x 3/2 (t 2 -ti) 



X~c(x). 



(5.30) 



Consequently, xz < x^, where x — xz is the solution of 
Xc(x) = 0, and x = xa is the solution of xc(x) = 0: 



xa = 



-1 2/3 



(ta - h) 



(5.31) 



VI. 



AN INITIAL VELOCITY AMPLITUDE 
CONSISTENT WITH THE CMB 
OBSERVATIONS. 



Recent observations of the power spectrum of the CMB 
radiation [l(| [HI, [HJ show that the maximal value of 
the temperature anisotropy, AT, is approximately 80fiK, 
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and is attained for perturbations with the wavenumber 
I ps 200, that corresponds to the angular size of nearly 1 
degre e fill . fl3l fll . Il5j ]. (The numbers in the graphs in 
Refs. [l3| - [15| imply that the "angular size" meant here 
is the radius rather than diameter of the perturbation.) 
Consequently, one might take 80fj,K as the upper limit of 
temperature variation and 1° as the upper limit on the 
angle subtended by the radius of the perturbation. How- 
ever, both these parameters far exceed the scale needed 
to model the formation of a galaxy cluster, and even more 
so for a single galaxy. 

For simplicity we take the k = Friedmann model as 
the background. The relation between the physical size 
at recombination, and the angle on the CMB sky was 
given in paper I. 

As shown in the table in paper I, the 1° angle on the 
CMB sky is more than 10 times the angle occupied by the 
mass that will later become an Abell cluster of galaxies, 
and about 250 times the angle occupied by the mass of 
a single galaxy. In order to get the right estimate for a 
galaxy cluster, one has to go down to the angular scales 
below 0.1°. In this region of the power spectrum observa- 
tional data are sparse and carry big error bars (see Refs. 
[l2l [Til]), but the amplitude is approx. 25/iK. At the 
angular scales corresponding to single galaxies, 0.004°, 
there is no good data. 

One way to get an upper bound on the fluctuations of 
velocity, is by assuming the observed temperature vari- 
ation is entirely due to local matter motion at recombi- 
nation causing a Doppler shift. Now, the cosmological 
redshift 



3V 
To 



So 



1 + Z 



(6.1) 



is due to the expansion of the space intervening between 
the emission point at recombination and the observation 
point at us today. A velocity fluctuation at recombina- 
tion — i.e. a fluctuation of the velocity of the emitting 
fluid relative to the average background flow — superim- 
poses an additional Doppler redshift 



AT 



v 



Av 

c 



(6.2) 



which is unchanged by the subsequent cosmological red- 
shift of both T and T + AT. Therefore AT/T ps 10~ 5 
implies Av/c ~ 10 . 

To properly estimate the magnitude of velocity fluc- 
tuations right after recombination, we use a paper by 
Dunsby [16| . His eq (84) relates the observed temper- 
ature fluctuations to a mode expansion of the velocity 
perturbations (3ll| . From this and the various definitions 
elsewhere in the paper, we find, for a single dominant 
mode of wavelength A and relative velocity Av, 



Av/c 



AT 16tt 3 ^ 
"T"3(5 t )3A3 



(6.3) 



Av/c = 0.001895. 



(6.4) 



and, assuming a wavelength scale containing an Abell 
cluster mass at recombination, ^5.87 kpc, 

AT _ 25fiK 
~T~ ~ 2.12K 

This is distinctly more than was estimated by assuming 
AT was entirely due to a Doppler shift, and the reason 
is that the effects (on the observed AT/T) of density, ve- 
locity and intrinsic temperature perturbations partially 
cancel each other. 

This must now be recalculated as a fluctuation in the 
expansion velocity, since our measure of velocity is b — 
Rj/M 1 / 3 , which is constant in space in the Friedmann 
limit. 

We use the k = FLRW dust model to estimate R t 
at recombination. For this model, E — and R(t, r) = 
rS(t), 

on 1/3 
9M (t-t B ) 2 



S 



and the past null cone is 



6 

' pnc ~ [W 
Thus 

R f t \pnc 2 



1/3 



(to ~ t B ) 1/3 - (t - t B ) 



1/3 



to — t B 



1/3 



- 1 



= 2 



t-t B 

= 2 [VT+~z - 1] . 

We can now write, for recombination, 

Ab/b = AR tt /R, t = (Av/c)/R, t 
AT 16n 3 Si 1 



- 1 



(6.5) 



(6.6) 



(6.7) 



T 3(5, t )?A 3 2[y/TT£ - 1] 
0.001895/(2[Vl + 1000 - 1]) 
9.896 x 10" 5 re 10~ 4 . 



(6.8) 



Returning to inhomogeneous models, what should the 
radius of the velocity perturbation be? The initial den- 
sity profile considered in paper I had an initial fluctu- 
ation of small mass which later accreted more mass to 
form the condensation. But with a velocity perturba- 
tion, accretion would be caused by a small initial inwards 
perturbation over the whole final mass. 

Nothing is known from observations about the possi- 
ble profiles of the initial velocity perturbation. Lacking 
this, we are free to choose the profile that will be easy to 
calculate with. 



VII. PROFILES & FITTING 

A. Geometrical Units 

For convenience of computation, geometrical units 
were employed in order to avoid extremely large numer- 
ical values. If a mass Mq is chosen as the geometrical 



12 



unit of mass, then the geometrical length and time units 
are L G = GM G /c 2 and T G = GM G /c 3 . It was conve- 
nient to define different units for different sized struc- 
tures. The units used for Abell clusters, and for voids 
are summarised in table [U 



Table Q] goes here 



mostly A2199 They give the mass at 0.2 Mpc and at 
1 Mpc for each cluster, estimated from fitting a model 
to the ROSAT X-ray luminosity data (see table HT|) . so 
this allows the determination of two parameters. We 
also attempted to control the maximum density and the 
compensation radius by trying profiles with a total of 4 
parameters, but these adjustments did not have a large 
effect overall. 



B. Background Model. 

We use the k — Friedmann model as our background 
or reference model, and all our compensation radii are 
calculated with respect to this model. It is chosen pri- 
marily for simplicity, and any other Friedmann model 
would do. 

For this model, with M as radial coordinate, the areal 
radius is 



2\ V3 



(7.1) 



where subscript b indicates the background rather than 
L-T functions or values, and the density is 



Table [TT] goes here 



We considered quite a variety of functional forms for 
our profiles. For a number of them we fixed the param- 
eters by fitting to the mass-radius data. A comparison 
of four of these fits with the UP (also fitted to the same 
data) is shown in figure [TJ and in table Hill Several runs 
were done with each of these. 



Table |TTT] goes here 



Figure [T] goes here 



Snp b = — 2 . 
The velocity (rate of change of R) is 

so that scaled velocity becomes 

, {Rt) b 



AfV3 



3t 



1/3 



(7.2) 



(7.3) 



(7.4) 



We could replace t with t—ts in all of the above, but we 
took tg = for our background model. 



C. Choice of Density Profile at t 2 = 1.4 ■ 10 10 yr 
1. An Abell Cluster 



2. A Void 

The main feature of a void is a low density region sur- 
rounded by walls of higher than background density. Our 
initial attempt to model a void was to choose a profile 
with p(M) increasing, but this was problematic because 
we did not have good control over how high the den- 
sity became before we reached the compensation radius. 
Thus we designed a profile with a density maximum at 
the wall, decreasing to background at large M. In this 
way no unreasonable density values are encountered even 
if the compensation radius is well beyond the density 
maximum. The profiles tried are summarised in table ITVl 
and figured 



Table HVl goes here 



The widely used 'universal profile' (UP) for the varia- 
tion of density with radius in condensed structures (e.g. 
[ill]), discussed in section [ill CI was found not suitable 
because of its divergent central density, and also it does 
not convert to p(M) in any nice analytic form[32|]. We 
therefore sought a profile that is similar in the region for 
which data are known, but had a finite — if large - 
central density, is expressed as p(M) and has an analytic 
R(M). 

To determine the parameters of a profile, we fitted it 
to the mass-radius data from [H for A2199 & A496, 



Figure [H goes here 



D. Choice of Velocity &: Density Profiles at 

ti = 10 5 yr 

For velocity fluctuations at recombination, we merely 
choose a harmonic wave of very low amplitude, with ei- 
ther 1 or 1.5 wavelengths across the diameter of interest. 
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The boundary of the region of interest is the radius that 
will become the compensation radius at time 2. See table 
IVland figure [3] for the details. 



Table [V] goes here 



Figure [3] goes here 



A few runs were done from initial density profiles, 
which are given in table IVT1 



Table IVT1 goes here 



VIII. PROGRAMS 

The foregoing procedure to solve for the L-T func- 
tions E(M) and ts{M) from given profiles R t ti(M) and 
p 2 (M), was coded in a set of Maple [2(| programs that 
also re-derived and verified the formulas, and fitted the 
profiles to the data. The programs to evolve the resulting 
L-T model and plot the results were written in MATLAB 
pll ] . Unlike the Maple programs, the latter did not need 
large changes from the Paper I versions, but were in any 
case refined and made more flexible, and in particular 
were adjusted to handle times where the big crunch has 
already happened in parts of the model. [33[ Also the den- 
sity to density solving program (Maple) was modified to 
handle both p\ > p 2 and p 2 > p%. 

As before, slightly different variables from those used 
in the foregoing algebra were more convenient for the 
numerics. In particular, for the velocity to density case, 
in place of x, a 2 , &i and we use y = xa 2 , fi = b 2 a 2 
and Zi ~ ti/a 2 3 / 2 . For the velocity to velocity case we 
use y — x/b 2 2 , 7 = b 2 /b 2 2 and Zi = tib 2 . Expressions 
for the following limits at the origin are also needed, 

2 { ' ' 

~ 6, 2 (0)t^(0)' 

where the upper sign is for the hyperbolic case and the 
lower one for the elliptic case; e = +1 is for expansion, 
and e = — 1 for collapse. Since x(0) comes from the solu- 
tion procedure, we have an expression for 0^(0), whether 
velocity or density profiles are given, but for &i(0) when 
the velocity profile is given, the limit must be calculated 



using the explicit velocity profile function. Alternatively, 
the velocity may be given in the form 6j(Af). 

For reconstructing the model evolution from the so- 
lution functions x and ts 1 the formulas in Paper I still 
suffice. 



IX. RUNS 

The runs that we did are listed in table lVIll which also 
lists the relevant figure numbers where appropriate, and 
the significance of each run. 



Table IVlTl goes here 



Run Figures go here 



X. RESULTS 

A run is defined by specifying one of the initial profiles, 
pi & R,ti, and one of the final profiles, p 2 & R.t2- The re- 
sults of a particular run are the two L-T functions E(M) 
and ts{M) — the local energy /geometry function and 
the local bang time function — that determine the L-T 
model. The model so defined may then be evolved over 
any time range, or its whole lifetime, and in particular, 
it can be verified that the given profiles are reproduced. 
In addition, of the 4 profiles pi, R.n, p 2 and R t t2, the 
two that were not used to define the run are also deter- 
mined by {E(M),tB(M)}, and are of interest, especially 
the initial one. 

• The primary Abell cluster forming run, V^i4/ofl6 
(fig. [5]) , used an improved present day density pro- 
file based on observational results, and also started 
from an initial velocity profile that was expanding 
less fast than the background at the centre, which 
is what one would expect to evolve into a conden- 
sation, see figs. [3] and O Though an Abell cluster 
model was generated in Paper I, the final profile was 
not so carefully chosen, and it was evolved from an 
initial density profile. However, the evolution of 
p, R and Rj is visually very similar to that given 
in paper I. The main differences are at t\. Sev- 
eral runs with different initial velocity profiles gave 
quite similar results — with only p\ and t b chang- 
ing significantly. 

• The primary void forming run, Ui3pf20 (figs. 
[5J [H [TO]), used an initial state that was expanding 
faster than the background at the centre, which is 
what one would expect to evolve into a low den- 
sity region. Using pf20 for the final density profile 
was particularly successful, as the required density 
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profile develops smoothly, and the density is still 
decreasing at every point even at the present time, 
so there is no sign of imminent shell crossings. 

It is evident from the plots of E & ts and the con- 
ditions for no shell crossings [6] , that there must be 
shell crossing at some time, but these are long be- 
fore recombination, and well after the present time. 
A significant feature is the relatively low central 
density at the initial time. The evolution of voids 
in an L-T model has been extensively discusse d by 
Sato and coworkers, see a summary in Ref. [23] ]. 
They, too, found that it necessarily leads to a shell 
crossing at some time. 

Our original attempt, using pl8, created a shell 
crossing (where the density diverges and then goes 
negative) so soon after t2, that the model evolu- 
tion program could only generate smooth density 
profiles up to 0.98t2- 

• The case of a pure density perturbation, ViOpilZ, 
was investigated by choosing a flat initial b\ evolv- 
ing to a cluster. This shows what initial density 
perturbation is needed at recombination to form 
the cluster by today — a central overdensity sur- 
rounded by an underdensity, of order 10~ 3 . 

Similarly, the case of a pure velocity perturbation, 
pi0pfl3, created by choosing a flat initial p\ evolv- 
ing to a cluster, shows what initial velocity per- 
turbation is needed — a general underexpansion of 
order 10 -4 , decreasing outwards. 

The latter is contrasted with the pure velocity per- 
turbation needed to create a void in pi0 / of20 - 
a general overexpansion of order 10~ 2 , decreasing 
outwards. (See fig. [5] for all 3 runs.) 

The rather large size of these perturbations tends 
to confirm that perturbations in both density and 
velocity are needed. 

• Our entire approach — specifying profiles for an L- 
T model at 2 different times — demonstrates that 
two L-T models with identical p(M) [and conse- 
quently R(M)] at a given time, or identical R :t {M) 
at a given time, can have quite different evolutions. 
This is shown by comparing runs pi0pfl3 & pi0pf20 
(both in fig [5]); also by runs VT4pfl6 & ViApiO; and 
in a small way by the set of 4 runs that all end 
with pfO. Thus, knowledge of only the density (or 
velocity) profile at one time does not determine a 
model's evolution at all. 

• To make clearer the effect of varying the initial ve- 
locity profiles, a set of runs was done with all 4 
initial V profiles and the homogeneous final p pro- 
file. We found ts has almost the same shape as the 
velocity perturbation &i; E is strongly influenced 
by bi, except that it starts from zero at M = 0, 
and returns close to zero at M = M c ; and the p\ 



perturbation is also strongly dictated by bi, but 
turns away from zero towards M = M c . 

• Run pi3pi2l (fig. [TT|) demonstrates the density pro- 
file inversion that [22j proved was not only possible 
but likely, and our method makes it easy to obtain 
such evolutions. 

This same run also shows a model can evolve from 
lower to higher density in one region and higher to 
lower in another; and that our method can handle 
such occurrences. 

• For the two times we have been considering — t\ 
at recombination and t% at the present day — the 
shape of the resultant E(M) in any given run is 
dominated by the final density profile, while both 
initial and final profiles have a noticeable effect on 
the resultant is(M). 

• We have found that, when choosing the initial ve- 
locity fluctuation with an appropriate magnitude, 
it is quite difficult to keep the initial density fluc- 
tuation (which is not chosen) small enough, and 
vice-versa. The main run of Paper I, which im- 
posed a final density profile of similar scale to pfl3, 
but less centrally concentrated, and a 3 x 1CP 5 den- 
sity fluctuation at ti, was quite successful in this 
regard, having a resulting velocity variation that 
was within 3 x 1CP 5 over much of its range, only 
reaching 8 x 10~ 4 in the outer regions that evolved 
into vacuum. This sensitivity to profile shape, and 
better way of choosing initial profiles could be in- 
vestigated further. 

• The compensation radius was close to a parabolic 
point in every run. This must in fact be the case 
to high approximation under the conditions we im- 
posed. Since we are using a k = Friedmann model 
as our reference or background model, then any L- 
T model which has the same M within the same 
R at the same t — ts must be locally the same 
(same R(t,M c )) as the background, and therefore 
parabolic. By definition, the compensation radius 
is where M & R are the same as in the background, 
but in general our method does not ensure t — ts 
matches there. At the beginning of the calculation, 
we do not know (t — ts) yet, it is one of the re- 
sults of the procedure. As explained is sec. IIVD1 
the individual values of the time coordinate have 
no meaning. 34] Only the differences, fa — t\), 
fa — ts), {t — tg), etc, appear in the calculations 
and are physical parameters. Nevertheless, we as- 
sume that our t^ equals fa — ts) in the background. 
In consequence, ts — in the background, while 
in general ts(M c ) 7^ in the perturbation. How- 
ever, in practice, our use of tiny deviations from 
the background value at recombination t\, ensures 
that the jump in tg can only be very slight. This 
is what ensures a close-to-parabolic model at M c . 
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As a result of the above, in each of the numeri- 
cal examples E(M) has the same sign in the whole 
range. However, the method discussed here and in 
Paper I can be freely applied to models in which 
E(M) changes sign at some M — see Paper I for 
an example. 



APPENDIX A: EXTENSIONS OF PAPER I — 
THE EVOLUTION BETWEEN TWO DENSITY 
PROFILES. 

1. Models with Both Larger and Smaller Densities 
at Later Times 



XI. CONCLUSIONS 

It is traditional to think of spacetimes as being 
specified by functions such as components of the 
metric and its derivative, or density and velocity, 
given on some initial surface. We have taken the 
new approach of specifying one function on an ini- 
tial surface and another on a final surface, and 
demonstrated that, in the case of the L-T model, 
an evolution from one to the other can be found. 
This approach is better suited to feeding the obser- 
vational data into the model functions. 



In particular, the spherically symmetric evolution 
from a given initial velocity (or density) profile to 
another given density (or velocity) profile, in any 
of the four possible combinations, can always be 
found. The solution is a determination of the arbi- 
trary functions that characterise the model. While 
there is no guarantee that the resulting model will 
be free of shell crossings, it is a simple matter to 
check the arbitrary functions once calculated, and 
our experience has shown that reasonable choices 
of the profiles keep any shell crossing well before 
the initial time or well after the final time. 



It is clear from the method of Paper I that, for any r 
value, whatever solution is found for p 2 < p\, its time 
reverse (t — ► t\ + 1% — t, ts — > t\ + t 2 — ts) will solve the 
case when the density values are interchanged. But could 
both p(ti,r a ) > p{t 2 ,r a ) and p{h,r b ) < p(t 2 ,r b ) occur at 
different r values r a & r b in the same model? It is obvious 
that a model with adjacent expanding and collapsing hy- 
perbolic regions must have severe shell crossings, unless 
they were separated by a neck [f|. But the latter puts 
the two hyperbolic regions on either side of a wormhole, 
so they don't really communicate. On the other hand it 
seems entirely possible that two worldlines in an elliptic 
region could display such behaviour, or even a hyperbolic 
region outside an elliptic region. 

In Paper I, the condition p 2 (M) < pi(M) was imposed 
merely to ensure R 2 > Ri, which allowed us to know 
there was always a bang in the past, and hence a well 
defined ts- However, we note that this condition did not 
exclude other equally unreasonable possibilities, such as 
an elliptic region outside a hyperbolic region (see Q). 

Thus, if we relax this requirement, and allow all p\ > 
and all p 2 > profiles, we merely have to note which of 
Ri & R 2 is larger before generating the "forwards" or 
"backwards" solution. 

We in any case have to check whether the conditions 
for no shell crossings Q are satisfied, or whether a regular 
maximum or minimum has been reached, and this merely 
adds one more thing to check — whether the bang and 
crunch functions are sufficiently continuous. 



The solution method has been programmed, and a 
variety of runs have demonstrated the practicality 
of the method. Models of a void and of an Abell 
cluster have been generated, and shown to have 
well behaved evolutions between the initial and fi- 
nal times. It is noteworthy that the void model was 
able to concentrate matter into a wall without the 
formation of any shell crossings up to the present 
time and for well after it. 



The effect on the solution of varying the initial 
and final profiles was effectively illustrated. Some 
other previously known features of the L-T model, 
such as the possibility that clumps could evolve into 
voids, were highlighted in the results section by pro- 
ducing examples. 



2. Including regions of zero density in the profiles 

a. Transient zeros in the 



Given the expression for the density, (|2.3[) it is clear 
that if M r = for a particular r, but E >r and (ts),r 
are not zero, the density will be zero there for all time. 
Is it possible that p = at isolated events, or on non- 
comoving worldsheets? 

Assuming M r ^ and < R < oo, the only way 
for this to happen would be for R r to diverge without 
changing sign, while Af r remains finite. We immediately 
see that divergent R_ r makes g rr divergent, which sug- 
gests bad coordinates at the very least. However, for all 
E values we may write @, HI 



R,r 

~R 



~E 



R 



(t-t B ) 



(A.l) 



R, 
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and it is clear that tB, r R,t is only divergent on the 
bang or crunch, while (i — tg)R t t is zero there, because 
R t ~ (t — is) -1 / 3 at early times. So one possibility 
is that t\ intersects a non-simultaneous bang while t 2 
does not intersect the crunch at the same r value, or vice 
versa. This of course means the constant time surface 
cannot be extended beyond this point, which is not re- 
ally satisfactory. The only other possibility is for E r or 
tB,r to diverge. [35]But this would give a comoving zero 
of density, and a coordinate transformation would make 
E r and/or tg r finite and M r zero. 

Thus we conclude that, apart from the case of the ini- 
tial or final surface intersecting the bang or crunch, the 
density along a particle worldline cannot be zero at one 
time and non-zero at another. Thus zeros of density have 
to be permanent and comoving. 



b. p(M) = at a single M value 

To look at the question of whether zero density at a 
single M value can be accommodated in the methods of 
Paper I, we first determine what kind of zero we might 
expect in p(M). Working in flat 3-d space, let 

p(R) w A\R-R z \ n , A,n const & > (A.2) 

near a zero in p. Then 



\M-M Z \ = 4tt / p{R')R' 2 dR' 



AttA 



\R-R : 



\n+3 



71 + 3 

Rl\R-RA n+1 



2R Z \R — R 2 
77 + 2 



\n+2 



n + 1 



so, to lowest order 



\M-M Z \ « AnA R l\ R - R *\ n+ \ 
n+1 

(n + l)\M - M z \\ 1/(n+1) 
A-kARI 



(A.3) 

(A.4) 
, (A.5) 



\R-R Z \ r 
which means 

and we notice 77/(71 + 1) < 1 always. 

If we use this as an approximation near a point zero in 
the full curved spacetime expression, we find 



R — Rt 



M 



G 



Mr 



np(M') 



■dM' 



(A.7) 



= ■^(n + l)\M-M z \ 1 ^ n+1 \ 



which is indeed well determined. We thus conclude that 
no modification of the Paper I method is needed in prin- 
ciple, though some extra coding might be needed if the 



integration had to be done numerically. (With all the 
profiles tried so far, Maple did symbolic integrations to 
get Rl) 



c. p — over an extended region 

We have already seen that zeros in the density have to 
be permanent and comoving, and this obviously applies 
to extended vacuum regions. In any case, spherically 
symmetric vacuum is Schwarzschild and must remain so. 

In order to allow the possibility that p(U, r) is zero over 
a finite range, the choice pi = p%(M) must be abandoned 
as all points in the zero density region will have the same 
mass, so M will be a degenerate coordinate there. 

Whatever alternative possibilities we choose must al- 
low us to identify the corresponding comoving coordinate 
points on each of the density profiles at t\ & t 2 - Thus 
R is not a suitable coordinate, as the two R(ti,r) are 
different. 

Thus we instead consider the following alternatives: 

1. Specify pi = p l {r) and M = M(r) 

This parametric version of choosing p(M), allows 
Pi(r) — and M(r) — const over some range of 
r, but actually is not sufficient, as we still have no 
idea how much Ri{r) increases over this range — 
how much space there is between the non- vacuum 
regions. In other words, in such a range, neither 
pi(r) nor M(r) provides a useable definition of r in 
terms of a physical quantity. 

2. Specify pi = pi(Ri) 
in which case 



M t (R) = (k/2) 



Rj Pt 



(R)dR. 



(A. 



Now where pi > 0, M(R) will be different for each 
R value, and therefore can be used to identify cor- 
responding points on the two profiles. We can then 
use i?2 or Ri (or something else) as the coordinate 
radius r. Let us say we choose r = R2. But for 
regions where pi = 0, M(R) will be degenerate. 
In this case, we take advantage of the fact that in 
vacuum there are no matter particles, and so we 
have extra freedom in choosing the geodesies that 
constitute our constant r paths. Using a linear in- 
terpolation between the R values at the two edges 
of the vacuum region at each of t\ & ti provides 
the obvious choice of corresponding points, i.e. at 
time t%. 



r = R 2 , 



(A.9) 



and at time t±, 
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r> | (Hi vac-fll vac min)(i?2 vac,max--R2,vac,min) 

i?2,vac,min + (R 1 ,:„-R 1 , mc , IaiTl ) Ul vacuum, 

MjT 1 (Mi(i?i)) in non- vacuum. 



The original procedure for extracting the functions 
E and is, which actually used R\ & R 2 , goes 
through with the only change that functions of 
r = M are now functions of r = R 2 . In particu- 
lar, in vacuum regions, each worldlinc has different 
pairs of R\ & R2 values, but the same M value. 



APPENDIX B: THE EVOLUTION OF R(t,M) VS 
THE EVOLUTION OF p(t,M). 

Let us take the relation between R and p in the most 
general case, without any assumptions, eq. (|2.10|) . 

This means that the value of R at any M depends on 
the values of p in the whole range [Mo , M] . Also the 
inverse relation, (|2.9[) . is nonlocal — to find p(t,M) we 
need to know R in an open neighbourhood of the value 
of M. 

Assume now P2 < p\ over the whole of [Mo, M\. Then 



R 3 (M) - R\{M) - (R 3 (t 2 ) - R Q 3 {h)) 



1 



1 



'Mo \P2(u) Pl(u) 

Hence, i?|(M) > R\{M) if 



du > 0. 



(B.l) 



Ro 3 (t 2 )-R Q 3 (t 1 ) + - 



M 



Mo 



1 



1 



p 2 (u) 



du > 0. 



This will hold if we assume, as was done in paper I, 
Ro(ti) — 0,z = 1,2 in addition to p 2 < P\- The assump- 
tion Mq = 0, i.e. the absence of a mass-point at R = 0, 
also made in paper I, is not needed for this purpose. 

The converse implication is simply not true: 
R(t 2 ,M) > R(h,M) for all M £ [M ,Mi] does not 
imply anything for the relation between p(t 2 ,M) and 
p(ti, M). This somewhat surprising fact is easy to under- 
stand on physical grounds. R(t 2 ,M) > R(t\,M) for all 
M £ [Mo, Mi] means that every shell of constant M = M 
has a larger radius at t 2 than it had at t\. However, the 
neighbouring shells may have moved closer to M at t% 
than they were at t\. If they did, then a local condensa- 
tion around M was created that may result in p(t 2 ,M) 
being larger than p(t\, M). This does not happen in the 
Friedmann limit, where local condensations are excluded 
by the symmetry assumptions. 

A sufficient condition for p 2 (M) < pi(M),M £ 
[M , Mi] is 



1 



1 



(Ri 



> 



'R2 3 ) 



for all M £ [M , Mi]. 



If R, M > for all M £ [M ,Mi] at both t x and t 2 (i.e. 
there are no shell-crossings in [Mo, Mi]), then this is 
equivalent to (i? 2 3 ) ,m > (Ri 3 ) ,m for all M £ [M , Mi]. 
Incidentally, this implies R 2 > Ri for all M £ [M ,Mi] 
ifJ2a(M ) >Ri{M ). 



APPENDIX C: CALCULATIONS FOR SECTION 

IV. 



1. Derivation of (|4.33p and (|4.34p 

In eq. (|4.31[) we have 

*x(0) 
$x (2/oa) = 



0. 



/« 2 6i 2 /2-l 



(Z2&172 + 1 

- (2/a 2 f 2 (t 2 - t x ) , 



\a 2 h 72 + I 



(C.l) 



«2 



3/2 



2&i 



V2 - a 2 x 
\fx\x {x). 



(C.2) 



+2/a 2 



We see that $x,k(0) = and that lim^ 
Now XxAx) = \a 2 5 ' 2 /{2 - a 2 xf' 2 + Ah/ih 2 + x) 3 > 
for all x £ [0,2/02]. Hence, Ax can have at most one 
zero and it will have one only if Xx(fi) < 0, i.e. if eq. 
(14733T) holds. 

With fulfilled, <5> x ,x < in a neighbourhood of 

x = 0, &x,x — at some x — x m i n , and &x,x —> +00 at 
x — > 2/02; the latter means that the tangent to $x(x) at 
x = 2/a 2 is vertical. Consequently, &x(x) itself is a de- 
creasing function for x £ [0, x m i n ), and is negative in this 
range, then it is increasing for x £ (x min , 2/a 2 ]. It can 
thus have a zero at any x > if and only if $^ (2/a 2 ) > 0, 
i.e. if (14341) is fulfilled. 



2. Derivation of (|4.40p . 

The relevant properties of $(7(0;) are 
$c(0) = 2tt > 0, 

a 2 6i 2 /2- 1 



_ ,„ , , 6iv2a2 

q>c(2/a 2 ) =7r H 5- arccos 

a 2 0i /2 + 1 

-(2/a 2 ) 3/2 (t 2 -h), 



026172 + 1 



(C.3) 
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< 0. 



a-2 



3/2 



2&i 



\/2 — a 2 x 



fa" 



+ o(*2-*l) 



(C.4) 



Since $c(0) > and <& c .x < for all x > 0, the function 
^ct^) may have (only one) zero anywhere in its range 
if and only if <I>c(2/a2) < 0, which translates into the 
opposite of (f434j) . i.e. eq. (f4~lQ]) . 



t 2 > t\. The proof comes out quite simply if we look at 
the function xx (x) defined by (|5 . 13[) as a function of the 
argument b 2 . At b 2 = b\ we have: 



Xx(x) =Xi{x) = -x 3 / 2 (t 2 - h) , (E.l) 



which is obviously negative for all values of x > 0. At 
b 2 - * 00 we have 



APPENDIX D: THE VELOCITY AT t 2 MUST BE 
SMALLER THAN AT ti IN THE HYPERBOLIC 
EVOLUTION. 

We now show that b 2 > b 2 2 must hold, as stated 
in subsection IV A 1| in consequence of t 2 > ii and 



R, 



!t=tl 



> 0, i.e. b\ > 0. (The assumption R it {t\) > is 



hidden already in eqs. - (EH); for R t < 0, (t - t B ) 

would have to be replaced by (tB — t) in all three places.) 
Suppose that b 2 2 > b\ 2 , so that x < b 2 applies. Then 



lim Xh(x) 



(D.l) 



and in addition \h (0) = (this second property does not 
depend on the sign of b 2 2 — b\ 2 ). It follows that a second 
zero of xh (x) will exist if there exists a subset of R\ on 
which xh is an increasing function. From Eq. (|5.6[) we 
see that Xh,x(0) = 0, lim^.^^ xh,x — —00. Hence, in 
order that xh(x) may have a zero at x > 0, the /ih(x), 
defined in (|5.6p , must be positive somewhere in the range 
x > 0. The derivative of hh{x) is 



fj,' H (x) 



4b 2 



46i 



(& 2 2 



and so there must exist values of x such that 



h (h 2 - 
This is equivalent to 

"(6i/6 2 ) 1/3 -l 



h (b 2 



> 0. 



(D.2) 



(D.3) 



x > 



6l V3 ^5/3 _ 6l 5/3j (D 4) 



However, this is a contradiction since, with 61 < 62 and 
x > 0, the left-hand side is negative, and the right-hand 
side is positive. 

Since b 2 2 > b\ 2 has thus led to a contradiction, the 
opposite must hold. This result is intuitively obvious 
(for dust, expansion must slow down with time), so the 
above is in fact just a consistency check. 



APPENDIX E: THE VELOCITY AT t 2 MUST BE 
SMALLER THAN AT t x IN THE ELLIPTIC 
EVOLUTION (FINAL STATE STILL 
EXPANDING). 



We will prove here that the b 2 in eq. (|5.13p must be 
smaller than b\ if, as was assumed earlier, b\ > and 



Xx (x) = X2(x) 



61 — X 



arccos 



hi 



x 3 ' 2 (h-h). 

(E.2) 



We find that X2(0) = 0, lim^oo X2i x ) — and the 
derivative is x' 2 i x ) = ~ 2b i Vx/(h 2 + x 2 ) 2 - (3/2)y/x(t 2 - 
t\) — obviously negative for all x > 0. Therefore X2(x) 
itself is negative for all x > 0. 



Now we calculate: 



dxxjx) 
db 2 



Ax 3 ' 2 



{b 2 



(E.3) 



which is negative for all values of b 2 at every value of 
x > 0. Thus, xx(x) is negative for all x > at b 2 = 61, 
negative for all x > at b 2 — > 00 and is a decreasing 
function of b 2 at every b 2 € (61, 00) for every x > 0. This 
means that xx{x) is negative at every value of x > 
for any 62 > 61, and so the equation xx{x) = has no 
solutions in (0, 00) when b 2 > b\. □ 
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More correctly, the elliptic, parabolic and hyperbolic con- 
ditions are E/M 2 ^ 3 < 0, = 0, and > respectively, as 
both E & M go to zero at an origin. 
One of the definitions of R m used in astronomy is: The 
size R m of a galaxy cluster is that radius, at which the 
average density within R m equals a specified multiple 
of the background density pb (e.g. 200 [||). With this 
definition, the procedure of determining \ii is different 
from the one presented in eqs. (|3.5[1 - (|3.6[) . 
A proper matching of first and second fundamental forms 
would have matched masses, energies and ages, and be 
a comoving surface M = M ma t c hed- If our 'compensa- 
tion' procedure were executed at each time, the resulting 
surface would not be properly matched, nor would it be 
comoving. 

Note that the relation between mass and density in the 
L-T model, eq. (12. 3p . is very different from the flat-space 
relation dM/d£ = 4irp£ 2 . We will carry over to relativity 
the equation p — p(M) resulting from the Newtonian re- 
lation with the hope that it is a good approximation at 
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low densities and for small radii from the centre. How- 
ever, a completely self-consistent approach would require 
re-interpretation of all the relevant astronomical obser- 
vations against the L-T model. 

This profile is meant for large £ only. However, all our 
numerical codes cover M = 0, so we added v to avoid 
modifying the codes. 

The L-T model is general enough to describe a situation 
where the final profile is also just a bump on a perfectly 
homogeneous background, but this picture does not fit 
well with the popular image of mass being accreted onto 
the initial condensation. The "bump on a smooth back- 
ground" at t2 would imply that the whole infinite back- 
ground adjusted its density to the central bump, which 
would imply propagating the density wave to an infinite 
distance in a finite time. 

In this way, we also avoid a numerical difficulty: the vol- 
ume that is numerically monitored must be finite. The 
"condensation surrounded by a rarefaction" model has 
the advantage that the region outside the rarefaction re- 
mains strictly Friedmannian for all time and it is enough 
to monitor it only at its edge. In the "bump on a smooth 
background" model, the velocity perturbation spreads 
out to infinity; only the mass-density at ti is Friedman- 
nian outside the bump (and even this homogeneity gets 
destroyed by evolution, both backward and forward in 
time — it is only momentary). 

The velocity used here is the emitting fluid velocity rela- 
tive to the normals to the constant temperature surfaces. 
We checked whether this profile would encounter a black 
hole horizon R — 2M. The mass within radius R is 
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47T q 

M = —p b 8R 3 s 



In 



l + R s 



1 + 



so while the function R — 2M could in principle have a 
root, it does not have one for reasonable astronomical 
parameter values. 

These will be useful in modelling the formation of black 
holes. 

Indeed, all physical quantities calculated are invariant 
under the transformation t — *• t+ const. 
If we allow M :r to diverge, then we need R' to diverge 
faster, and the argument is pretty much the same. 
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TABLES 

TABLE I: Geometrical Units. The sets of geometrical units used for models of Abell clusters and for of voids. 



Object 


M G 


Lg 


T G 


Abell clusters 


10 ib M 


47.84 parsecs 


156.0 years 


Voids 


55 x 10 ib Mq 


2631 parsecs 


8582 years 



TABLE II: Abell Cluster Data. Data from [H], giving the observed masses within 0.2 Mpc and 1 Mpc deduced from 
X-ray observations, which is used to fit the profiles in table ITTTl to observations. 



Abell Cluster 


M(0.2 Mpc) 


M{1 Mpc) 


A 2199 


(0.65 ±0.11) x 10 14 M e 


(2.9 ±0.3) x 10 14 Mq 


A 496 


(0.47 ±0.10) x 10 14 M e 


(3.1 ±0.3) x 10 14 M e 



TABLE III: Abell cluster profiles. The equations and properties of the profiles shown in fig [T] that were used to 
model the present density distribution in Abell clusters. The parameters given are those that fit the A2199 data Q3 . All 
non-dimensionless values are in geometric units — see table [I] 



Profile 


p/pb 


Parameters 


Rc 


M c 


pfUP 


P _ 4 


8 = 77440, 
R« = 3457 Lg 


324600 L G 


0.9442 M G 


pfl3 




B 2 = 7774000, 

\i% = 0.05304 Vic , 

(J 2 = 2/l 2 , 
^2=5 


374500 L G 


1.449 M G 


pfl5 


p B^- M l<>1 


B 2 = 507500, 

m = 0.04474 M a , 

(72 = 2/i 2 , 
V2 — 5 


312400 L G 


0.8416 M G 


pfl6 


P _ B 2 
Pb l+e M /M 


B 2 = 103400, 

A*2 = 0.04577 VMg 


286600 L G 


0.6500 M G 


pfl7 


P _ B 2 
Pb l-)- e VM / l i 2 


B 2 = 498500, 
p 2 = 0.07144 M G 


349800 L G 


1.182 M G 



TABLE IV: Void Profiles. The equations and properties of the profiles shown in fig [2] that were used to model the 
present density distribution in voids. The parameters were chosen to give a low density region well below the background 
value (see [p/pb]o) within 60 Mpc radius. All non-dimensionless values are in geometric units — see tableU 



Profile 


P/pb 


Parameters 


[p/Pb]o 


Rc 


M c 


pfl8 


= A 2 + B 2 e M/ ^ 

Pb 


A 2 = 0.01, 
B 2 = 0.01, 
p 2 = 0.01429 M G , 


0.01 


22804 L G 


0.9902 M G 


pf20 


p _ B^p^ + M'-' 


B 2 = 0.3873, 
C 2 = 1, 

p 2 = 0.2475 M G , 


0.075 


35850 L G 


3.846 M G 


Pb ~ C 2 p'i + (M-p 2 ) 2 



TABLE V: Initial Velocity Profiles. The equations and properties of the profiles shown in fig [3] that were used to 
model the initial velocity fluctuations at recombination. The values of M c were determined by the final density profile being 
used in each run. All values are dimensionless. 



Profile 


R,t/(R,t 


)b 


Parameters 


ViO 




1 




(R.t)b 




Vi3 


R,t 


1 + Ai (l + cos(^)) 


Ai = 5 x 10" 5 , 
M c from t 2 data, 


(R,t)b 



ViA 


T «^ = l-A 1 (l + cos( 


kM \ \ 
M a )J 


Ai = 5 x 10 -5 , 
M c from t2 data, 


Vi5 


(nX- 1 + AlCOS 


V 2Af c , 


) 


Ai = 1 x 10" 4 , 
M c from t2 data, 


Vi6 




' 3ttM 
^ 2M C / 


) 


Ai = 1 x 10 -4 , 
M c from t2 data, 



TABLE VI: Initial and Final Density Profiles. The equations and properties of the profiles that were used to model 
the initial density fluctuations at recombination and the final structure. 



Pro- 
file 


P/Pb 


Parameters 


Description 


pfO 


P/Pb = 1 




Uniform background density 
at time 2. 


piO 


P/Pb = 1 




Uniform background density 
at time 1 . 


pi3 


P/Pb l-Ai(l-2M/ni) 


Ai =0.1, 
pi = 1. 


Central overdensity plus sur- 
rounding underdensity of or- 
der 10 _1 in region of mass 
1O 15 M . Used for the initial 
profile in run pi3pf21. 


pi4 


/- 

1 


IPb = 

' 1 + (l + cos(i2^)), 


Ai = 1.5 x 10 -5 , 


Overdensity of order 3 x 
10 -5 in region of mass 
1O 13 M , with exactly back- 
ground density outside it. 


pf21 


-/„ _ 1 Pbl 
P/Pb i+A 2 (l-2M/ M2 ) p b2 


A 2 =0.1, 

P2 = 1. 


Central underdensity plus 
surrounding overdensity of 
order 10 _1 in region of 
mass 1O 15 M0, with a 'back- 
ground' density equal to 
that at time t\. Used for the 
final profile in run pi3pf21. 



TABLE VII: Modelling Runs. 



Initial 
profile 


Final 
profile 


Description 


Function (Fig 
no) list 


What it shows 


pi4 


pfl7 


small localised density per- 
turbation — > Abell cluster 


E(M) d4J), t B (M) ((4J), 

R,ti( M ) ©, 
P (t,M) ©, 


Classic evolution of a small 
initial fluctuation with 100th 
final mass into a conden- 
sation with realistic density 
profile at £2- 


ViA 


pfl6 


Low central initial expan- 
sion rate — > Abell cluster 


E{M) ©, t B {M) ©, 
pi(M) © 


Evolution of low expansion 
region, to a realistic density 
profile at t^- Graphs fairly 
similar to Paper I model. 


ViO 


pfl3 


Uniform initial expansion 
rate — > Abell cluster 


£(M)(6] left), t B (M)® 
left), pi(M)dl left), 
p(i,M)P left), 


Case of a pure density per- 
turbation — how big a per- 
turbation is needed for a 
present day structure? 


piO 


pfl3 


Uniform initial density — > 
Abell cluster 


B(M)P centre), 
t B (M)® centre), 
-R,ti(M)(H centre), 
p(t, M)(6] centre), 


Case of a pure velocity per- 
turbation — how big a per- 
turbation is needed for a 
present day structure? 


piO 


pf20 


Uniform initial density — > 
Void 


£(M)P right), 
t B (M)® right), 

R,ti(M)B "ght), 
p(t,M)(H right), 


Together with pi0pfl3, shows 
two models with the same 
density profile at one time 
can have totally different 
evolutions. 
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Vi3 


pfl8 


High central initial expan- 
sion rate — > Void 


(None) 


Not so good example of evo- 
lution to a void — it is so 
close to a shell crossing at 
t2 that the evolution program 
fails on the last time step. 


Vi3 


pf20 


High central initial expan- 
sion rate — + Void 


E(M) J7|), t B {M) (8J, 
pi(M) ©, p(t,M) HO]), 


Example of evolution to a 
void with good behaviour at 
£2 and well past it. The den- 
sity perturbation at t\ is too 
large by 4 orders. 


Vi3 


pfO 


High central initial expan- 
sion rate — > Uniform density 


(None) 


These 4 show very clearly the 
effect on E(M) & t B (M) of 
varying initial density 
or initial velocity. Since the 
final density is uniform, the 
deviation of E(M) & 
t B (M) from their FLRW 
forms is entirely due to non- 
uniform initial profiles. 
They also show a given 
density profile at one mo- 
ment can have many different 
evolutions. 


Ui4 


pfO 


Low central initial expan- 
sion rate — > Uniform density 


Vi5 


pfO 


High central & low outer ini- 
tial expansion rate — > Uni- 
form density 


Vi6 


pfO 


Low central & high outer ini- 
tial expansion rate — > Uni- 
form density 


pi3 


pf21 


Expanding underdensity — > 
collapsing overdensity 


E(M) t B {M) 1UI, 

pi(m) dnj, p 2 {M) cehi, 


An example of p2 < pi 
and within the same model 
pi < p2- It also provides 
an example of density profile 
inversion. 
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FIG. 1: Comparison of density profiles used for Abell clusters. The curves are identified by the numbers 
given in table [TTTI Although each curve is matched to the same pair of (R, M) points (see table[n]l, this only 
ensures the average density within each of those two radii is the same for all profiles, not the actual density. 
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FIG. 2: Comparison of density profiles used for voids. The curves are identified by the numbers given in 
table |IV] 
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FIG. 3: Comparison of initial velocity profiles at recombination. The curves are identified by the numbers 
given in table [V] 




FIG. 4: E(M), ts(M), the R } tl(M) fluctuation, and the p(t, M) evolution for run pi4pfl7. Notice the effect of the central bump in 
pl(M) shows up in the tg and R t ti curves. Although dtg/dM > indicates there will be a shell crossing, the small magnitude of tg 
ensures it will occur long before ti , when the model first becomes valid. 




' 0.1 0.2 0.3 0.4 0.5 0.6 0.7 

M 

0.7 1 1 1 1 i — 1 1 1 




0.1 0.2 0.3 0.4 0.5 0.6 0.7 

M 

FIG. 5: E(M), ts(M) and the pi(M) fluctuation for run Vi4pfl6, in which a region of initially low 
expansion rate evolves into an Abell cluster. In this and other figures, "H" , "EX" , "EC" , "Pe" and "Ph" 
indicate regions that are hyperbolic, elliptic and still expanding at time <2, elliptic and recollapsing at time 
t2, elliptic but within the range for a scries expansion about the parabolic model, and hyperbolic but within 
the range for a series expansion about the parabolic model. 
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FIG. 6: E(M), t B {M) the pi(M) or R tt i(M) fluctuation, and the p(t,M) evolution (top to bottom) for runs Vi0pfl3, pi0pfl3 and 
pi0pf20 (left to right). The left and middle columns show the pure density and pure velocity fluctuations needed at recombination 
to create an Abell cluster today. The middle and right columns show the difference in the pure velocity fluctuations needed at 
recombination to create an Abell cluster and a void today. 
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1 



0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 

M 

FIG. 7: E(M) for run Vi3pf20, which successfully creates a present day void. Though shell crossings arc 
inevitable, they don't occur for a long time after the present. 
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M 

FIG. 8: t B (M) for the void run Vi3pf20. 
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M 

FIG. 9: pi(M) fluctuation for the void run Vi3pf20. The amplitude of this fluctuation is too large, so the 
initial velocity profiles need finc-tuninng. 








FIG. 10: p(t, M) evolution for the void run Vi3pf20. Note that p2(M) is the same profile as p2(R) shown 
in figure [2] though it looks quite different. The reason is that there is very little mass in the void interior, 
so a large increase in distance corresponds to a small increase in mass. 
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FIG. 11: E(M), ts(Af) and the pi(M) and p2{M) fluctuations for run pi3pf21. The densities pi & p2 at times t\ & ti are of 
comparable magnitude (which makes p2 around 2 X 10 10 times the parabolic background at t\). In particular, pi > p2 at the centre, 
M = 0, and p2 > pi at the edge, M = 1. 



